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Controlling transitions in a Duffing oscillator by sweeping parameters in time. 
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We consider a high-Q Duffing oscillator in a weakly non- linear regime with the driving frequency 
cr varying in time between g% and <jf at a characteristic rate r. We found that the frequency sweep 
can cause controlled transitions between two stable states of the system. Moreover, these transitions 
are accomplished via a transient that lingers for a long time around the third, unstable fixed point of 
saddle type. We propose a simple explanation for this phenomenon and find the transient life-time 
to scale as —(In \r — r c |)/A r where r c is the critical rate necessary to induce a transition and A r is 
the repulsive eigenvalue of the saddle. Experimental implications are mentioned. 

PACS numbers: 05.45.-a 



Even the simplest of non-linear dynamical systems, 
such as those described by second order ODEs are noto- 
rious for exhibiting rich phenomenology One of the 
features of such systems is multi-stability, which strictly 
exists when all parameters are time-independent. In the 
case of adiabatically varying parameters, a system ini- 
tially situated at one of the quasi-fixed points will re- 
main close to it. What happens when the parameters 
vary faster then the time scales determined by eigenval- 
ues around the fixed points (FP) will be explored here on 
the case of a damped, driven Duffing oscillator. When 
the variation of parameters is sufficiently rapid, transi- 
tions occur, such that after the rapid variation is finished, 
the system finds itself at the FP different from the one 
at which it was initially situated. 

Statement of the Problem. We will consider an oscilla- 
tor obeying the following equation of motion: 



2/ cos 0(f) 



(l + ax 2 ) (1) 

7(f) 



where ujq is the angular frequency of infinitesimal vibra- 
tions, A is the damping coefficient, which for a quality 
factor Q is given by A = a is a non-linear coefficient, 
/ is the driving strength, and 7 is the driving frequency. 
Upon non-dimensionalizing and re-scaling we obtain the 
following equation which has the form of a perturbed 
simple harmonic oscillator: 



x + x = e (2F cos <p(t) 
<j> = l + ea(i) 



x — X 



(2) 



For the case of time-independent parameters F and a 
the amplitude response has the well-known frequency- 
pulled form 0, 0. For F > F cr there is a region of 
tri-valuedness with two stable branches and one unsta- 
ble (middle) branch (Fig. QJ. Such response curves have 
recently been measured for NEMS 4], indicating their 
Duffing-like behavior. For a constant a but different F 



there is also a response function with tri-valuedness for 
a > a cr . To each case there corresponds also a tri- valued 
phase response (not shown). 
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FIG. 1: Left: constant F. Right: constant a. 

We would like to explore what happens if the driving fre- 
quency, starting from the single- valued region (at c^) is 
rapidly varied in time, into the tri- valued region (ending 
at a/), at a constant F. We will also explore the phe- 
nomenology resulting from varying F at a constant a in 
a similar manner - from a single- valued regime into the 
tri- valued regime. The variation will take place via a 
function with a step - either smooth (such as tanh) or 
piece- wise linear - Fig. |3 In this article we will explore 
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FIG. 2: cr(t). 

the situation when F/F cr ~ 0(1), for instance F — 2F cr . 
When F ~ F cr , the size of the hysteresis ((72 — <ti) is also 
of order 1 and the perturbation method is especially sim- 
ple. The more complicated case of F ^> F cr (but small 
enough that r.h.s. of Eq. J2J) is still a perturbation over 
harmonic oscillator equation) may be considered in a fu- 
ture work. 

Phenomenology. We discover the following: upon 
sweeping a from lower values into the hysteresis, depend- 
ing on other parameters, the solution of Eq. (J2J) may jump 
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unto the lower branch. The transitions have a peculiar 
feature of having lifetimes much longer then the slow time 
scale of Eq. (0) - the damping time-scale of order e. 
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Transient with lifetime » run-down time 



FIG. 3: Example of a transition. 

Now, there are four relevant "control knobs": F, 07 
and r defined here as 1/ (sweeping time). Consider a set 
of imaginary experiments, each experiment performed for 
a different sweep rate r, with all three other parameters 
fixed. As r gets larger, depending on the value of other 
parameters there may be a critical sweep rate r cr beyond 
which transitions will be induced. Moreover for r ap- 
proaching very close to r cr the life-time of transitions, 
r, will grow. We will sample these lifetimes and r — df 
configuration space numerically and explain the finding 
theoretically below. However one more observation, de- 
picted in Fig. El is in order. For r close to r cr , not only 
do the life-times, r of transients tend to grow, but also 
the amplitudes and phases of these long transients tend 
to approach that of the middle (unstable) branch of the 
static Duffing oscillator. 
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FIG. 4: Solutions kiss the unstable branch - envelope of the 
full numerical solution of Eq. JJJ along with quasi fixed points 
(solutions to Eq. ©-0 below with lhs set to 0) are shown. 
The function a(t) is displayed as a dotted curve. 

We learn that for r close r cr the solution moves unto 
the unstable branch, lives there for a time period r and 
then either moves unto the top branch if r < r cr or per- 
forms the transition unto the bottom branch if r > r cr . 



The jump unto either the top or the bottom branch takes 
place long after reaching the static conditions (see Fig.[3J). 
A numerical experiment performed at particular param- 
eter values (see Fig. 0) demonstrates the typical situa- 
tion: plotted on the semi- log plot, the life-time versus 
\r — r c \ nearly follows a straight line. Thus, we learn that 
r oc — In \r — r c \. 
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FIG. 5: Transient lifetime VerSUS \V Tcr 

I for F « 2.8F cr , 

Gi — -2, Gf — 0.5((Jl + cr 2 ) 

Another numerical experiment was performed to measure 
r c versus 07 for F ~ 2.8F cr and <ii = —2 (see Fig. EJ) - 




FIG. 6: F = 2.8F cr , Gi — —2. There is a singularity at 
G l f. The size of the hysteresis was calculated with the first- 
order theory and therefore it is an over- estimation - the point 
where r c goes to zero (i.e. at infinitessimal sweep rates) is 
the exact right edge of the hysteretic region. The existence of 
singularity at g l f is a genuine effect, since the first-order theory 
predicts the lower end of the hysteresis more accurately. 

Two interesting features are immediately observed: (i) 
the critical rate necessary to induce a transition is not 
simply determined by the small parameter e in the 
Eq. J2J); (ii) when sweeping not deeply enough into the 
tri-valued region no transitions can be induced for any 
value of the sweep rate. Increasing <ii moves the singular- 
ity of r((jf) to larger values and varying F does not signif- 
icantly affect the shape of the r(cr/) curve. It is somewhat 
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surprising to see singularities in this rather simple prob- 
lem, especially before the onset chaotic regime. Transi- 
tions may occur by sweeping a down; also by sweeping F 
and holding a constant. Transition phenomenon persists 
for F large enough that the hysteresis is large (see note 
in Discussion below). 

Theory for small Aa case. The regime of F/F cr ~ 
O(l) guarantees that (02 — <ti) ~ 0(1), thus if is chosen 
sufficiently close to the hysteretic region, the jump in fre- 
quency during the sweep, Act is also ~ 0(1). This paves 
way for a simple perturbation method. For example, we 
write the "multiple-scales" perturbative solution to Eq. 
(2) as X = X (t,T, ...) + eXi(t,T, ...) + 0(e 2 ) where the 
slow time scale T = et, and in general = e n t. Plug- 
ging this into the equation and collecting terms of ap- 
propriate orders of e teaches us that xo = A(T)e lt + c.c, 
i.e. the solution is essentially a slowly modulated har- 
monic oscillator, with the modulation function satisfying 
the following Amplitude Equation (AE): 

Fe iS(T) _^_ iA _ 3 | A |2 A = Q 

where $ = ea(t(T)). Such AE holds for any sweep 
rate as long as Act ~ 0(1). The AE, broken into real, 
x = Re[A], and imaginary, y = Im[A], parts are: 

^ = -Ix + aWy-ltf + tfy (3) 

% = ~<T)x - \y + I (x 2 + y 2 ) x - | (4) 

These AE are well known and appear in similar forms in 
literature 0, . For a certain range of parameters these 
equations give rise to a two-basin dynamics with a sta- 
ble fixed point (FP) inside each basin. The basins are 
divided by a separatrix which happens to be the stable 
manifold of the unstable FP of saddle type. 
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FIG. 7: AE dynamics - an example. Typical solutions are 
shown as dashed curves. Black dots indicate fixed points. 

Such basins have recently have been mapped for a certain 
type of NEMS 6]. Next, consider a set of thought experi- 
ments, each sweeping a at a different sweep rate r. Imme- 
diately after the sweep, the system will find itself some- 
where in the two-basin space corresponding to conditions 



at a/, call it a point (xf(r),yf(r)). If (x f (r) , y f (r)) lies 
in the basin which has evolved from the basin in which 
the system began before the sweep then there will be 
no transition. If (x f (r) , y f (r)) lies in the opposite basin, 
then that corresponds to a transition! For very low r, 
during the sweep the system will follow closely to the 
quasi-FP and there will be no transition. In the opposite 
extreme - infinite r, at the end of the sweep the system 
will not have moved at all, due to continuity of a dynam- 
ical system. This endpoint of {(xj(oo), yf (oo))} may lie 
in either basin depending on the erf. 




FIG. 8: Real examples at F = 2.8F cr , en = -2. The locus 
{xf,Vf} is shown by dotted curves. 

In the numerical experiments that we considered the first 
point of the {(x f (r) , y f (r))} curve to cross the separatrix 
happens to be this tail 8] at r = oo. This explains the 
singularity mentioned earlier 9] (see Fig.EJ). 

Calculation of lifetimes. The point (x f (r) , y f (r)) 
serves as an initial condition for subsequent evolution 
at fixed Of. The (x f (r) , y j (r)) close to the separatrix 
(i.e. for r close to r c ) will flow towards the saddle and 
linger around it for a while. Because this lingering will 
take place close to the saddle, the linearized dynam- 
ics around the saddle should be a good approximation: 
r(T) = 5lv r e XrT + Rv a e XaT ', where, for example, v r 
and A r are a repulsive eigenvector and eigenvalue, SI is 
the distance of (x f (r) , y f (r)) away from the separatrix 
along v r and R is the characteristic radius of lineariza- 
tion around the saddle. The times at which the system 
crosses this circular boundary is given by R = yfr • r. 
The first time, Ti n , is of the order of SI. The second 
time is T out « j- In (J|) (neglecting effects of v r ). The 
lingering time is r = T out - T in w T out . So, 

1 . fSl\ 1 ( f r \r-r c \ \ ^ 

Thus we capture the logarithmic dependence of the tran- 
sient time versus \r — r c \. We found that calculating the 
eigenvalues A r using the second order theory as described 
in H systematically lowers the discrepancy in the slope 
of r vs. \r — r cr \ computed from exact Duffing and the- 
oretical predictions, Fig. which leads us to conclude 
that the errors are due to inexactness of the first order 
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relative error 
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FIG. 9: Dotted lines: Fractional errors between the slope of 
r vs. In \r — r cr \ from exact Duffing and l/(eA r ) where A r 
are computed using the first order perturbation theory in e 
(Eqns. (2)- (4)). Dashed lines: same comparison but with A r 
computed to second order. In both cases, F = 1.5, F = 2, 
F = 2.5 from left to right. Solid lines: the comparison of 
l/(eA r ) with the same slope from the AE, not full Duffing 
(only for F — 2 and F = 2.5) - remaining errors may be due 
to inexactness of linearized approximation and omission of v r . 

AE, not due to the incorrectness of the explanation of 
the cause of the transition phenomenon. 

Experimental Significance. One can propose to use the 
transition phenomenon with its long transients to po- 
sition Dumng-like systems unto the unstable (middle) 
branch (see Fig. [IJ - the desire to do this has been ex- 
pressed by workers in the NEMS community |7| . The first 
question is whether a necessary r c is attainable. We see 
from Fig. and related discussion that for a vast range 
of parameters r c is less then 0.1. Recall that in this 
paper r is defined simply as l/(sweep time) = 1/AT. 
The more experimentally relevant quantity is Aa /AT. 
In the present paper we are concerned with hysteresis 
widths of order 1, so Aa in question is ~ 10 or less. 
Hence Aa/AT necessary to create a transition is, for a 
vast range of parameters, less then 1 resonant widths per 
run-down time (but close to the high-end of the hysteresis 
this figure falls rapidly - see Fig.[5J), which in conventional 
units corresponds to the sweep rate of (coo/Q) 2 Hz/ sec. 
One can also use this as a guide to prevent unwanted 



transitions in experiments in which a depends on time. 
The second question is how small must \r — r c \ be to in- 
duce a transient of time r. From Eq. (??) we see that 
\r-r c \ « ^(dr/dl)xexp[-T\ r (<Tf,F)] sec -1 . For any F, 
A r reaches maximum approximately in the middle of the 
hysteresis where X r (a f = « .64(F - F cr ) / F cr . 

So the smaller the F, the easier it is to attain a longer 
transient. The quantity dr/dl, at the point of crossing 
the separatrix diverges at 07 = a 1 ^ and becomes small 
(< 1) for (jf close to a 2 (this explains why the life-time 
r is very sensitive to \r — r c \ in this region (see Eq. JSJ)). 
The functional form of dr/dl versus parameters is not 
fundamental - it depends on the form of the ramping 
function, for example, and can be computed from the set 
{x f (r),y f (r)}. 

Discussion. Our theoretical hypothesis claims that 
the point at which the system finds itself at the end 
of the ramp, (xf,yf) determines whether there will 
be a transition or not. However one may be inquire 
whether this set of points does not consistently cross 
the separatrix in some special way, for example, always 
tangentially or always perpendicular to the separatrix. 
If it does, there must be something happening during 
the ramp that situates these final points in this special 
way. Then the theory based on just the final position 
although is true, would be incomplete. This 
issue may be addressed in a future work by analyzing 
the effect of a generic perturbation of either the ramping 
function or the system on This might also 

pave the way to understanding the generality of the 
phenomenon - whether it holds in a large class of 
two-basin models and ramping functions. Also, we 
would like to explore the phenomenology for large Acr 
case when the AE (0-® do not hold, yet the system is 
still in a weakly non-linear regime. 
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